Trabajo Practico N° 1

Elaborado por Ramiro Duarte y Gonzalo Machin

se puede acceder a este trabajo practico tambien en el siguiente Link

Analisis de viajes en omnibus en la ciudad de Montevideo

Librerias

library(tidyverse)

## -- Attaching packages --------------------------------------- tidyverse 1.3.2 --
## v ggplot2 3.3.6     v purrr   0.3.4
## v tibble  3.1.8     v dplyr   1.0.8
## v tidyr   1.2.0     v stringr 1.4.0
## v readr   2.1.2     v forcats 0.5.1
## -- Conflicts ------------------------------------------ tidyverse_conflicts() --
## x dplyr::filter() masks stats::filter()
## x dplyr::lag()    masks stats::lag()

library(lubridate)

## 
## Attaching package: 'lubridate'
## 
## The following objects are masked from 'package:base':
## 
##     date, intersect, setdiff, union

library(sf)

## Linking to GEOS 3.10.2, GDAL 3.4.1, PROJ 7.2.1; sf_use_s2() is TRUE

library(mapview)
library(maptools)

## Loading required package: sp
## Checking rgeos availability: TRUE
## Please note that 'maptools' will be retired by the end of 2023,
## plan transition at your earliest convenience;
## some functionality will be moved to 'sp'.

library(ggrepel)
library(patchwork)
library(rayshader)
library(osmdata)

## Data (c) OpenStreetMap contributors, ODbL 1.0. https://www.openstreetmap.org/copyright

library(ggmap)

## Google's Terms of Service: https://cloud.google.com/maps-platform/terms/.
## Please cite ggmap if you use it! See citation("ggmap") for details.

Dataset de Viajes

Primero descargamos datos de los viajes realizados en el sistema de transporte metropolitano (STM) en Montevideo, como el archivo es grande lo descargamos y descomprimimos directamente de la pagina del catalogo de datos abiertos

temp <- tempfile()
download.file("https://catalogodatos.gub.uy/dataset/b1b22d81-9333-4a1b-8254-589268a698bf/resource/aa84ab90-7934-49ae-a330-57403f7e4e2e/download/viajes_stm_072022.zip",temp)
datos_julio <- read.csv(unz(temp, "viajes_stm_072022.csv"))
unlink(temp)

SHPs a Utilizar

Cargamos Shapefiles con recorridos de los omnibus, paradas y Zonas censales para contextualizar geograficamente los datos previos dado que no se cuenta en los mismos con un x e y

zonas_censales <- st_read("data/censo/ine_seg_11.shp")

## Reading layer `ine_seg_11' from data source 
##   `E:\Gonzalo Machin\Maestria\Materias\CDC 2\tp1_ciencia_de_datos_2\data\censo\ine_seg_11.shp' 
##   using driver `ESRI Shapefile'
## Simple feature collection with 4313 features and 13 fields
## Geometry type: MULTIPOLYGON
## Dimension:     XY
## Bounding box:  xmin: 366582.2 ymin: 6127919 xmax: 858252.1 ymax: 6671738
## CRS:           NA

paradas <- st_read("data/paradas_shp/v_uptu_paradas.shp", options = "ENCODING=UTF-8")

## options:        ENCODING=UTF-8 
## Reading layer `v_uptu_paradas' from data source 
##   `E:\Gonzalo Machin\Maestria\Materias\CDC 2\tp1_ciencia_de_datos_2\data\paradas_shp\v_uptu_paradas.shp' 
##   using driver `ESRI Shapefile'
## Simple feature collection with 38764 features and 10 fields
## Geometry type: POINT
## Dimension:     XY
## Bounding box:  xmin: 554395.5 ymin: 6134611 xmax: 591653.9 ymax: 6158038
## Projected CRS: WGS 84 / UTM zone 21S

recorridos <- st_read("data/recorridos/v_uptu_lsv.shp")

## Reading layer `v_uptu_lsv' from data source 
##   `E:\Gonzalo Machin\Maestria\Materias\CDC 2\tp1_ciencia_de_datos_2\data\recorridos\v_uptu_lsv.shp' 
##   using driver `ESRI Shapefile'
## replacing null geometries with empty geometries
## Simple feature collection with 652 features and 8 fields (with 2 geometries empty)
## Geometry type: LINESTRING
## Dimension:     XY
## Bounding box:  xmin: 554322.5 ymin: 6134518 xmax: 591684.6 ymax: 6158125
## Projected CRS: WGS 84 / UTM zone 21S

recorridos_no_maximales <- st_read("data/recorridos_no_maximales/uptu_variante_no_maximal.shp")

## Reading layer `uptu_variante_no_maximal' from data source 
##   `E:\Gonzalo Machin\Maestria\Materias\CDC 2\tp1_ciencia_de_datos_2\data\recorridos_no_maximales\uptu_variante_no_maximal.shp' 
##   using driver `ESRI Shapefile'
## Simple feature collection with 710 features and 13 fields
## Geometry type: LINESTRING
## Dimension:     XY
## Bounding box:  xmin: 555656.7 ymin: 6134558 xmax: 591668.1 ymax: 6158125
## Projected CRS: WGS 84 / UTM zone 21S

Correcion de Datos

Corregimos crs y cambiamos todos para trabajar unicamente en CRS = 4326

zonas_censales <- zonas_censales %>%
  st_set_crs(32721) %>%
  st_transform(4326)

paradas <- paradas %>%
  st_transform(4326)

recorridos <- recorridos %>%
  st_transform(4326)

recorridos_no_maximales <- recorridos_no_maximales %>%
  st_transform(4326)

Filtramos los datos de las zonas censales para unicamente trabajar con los datos del departamento de Montevideo

zonas_censales <- zonas_censales %>%
  filter(NOMBDEPTO == "MONTEVIDEO")

Visualizacion de Datos Espaciales

Observamos los datos geograficos al momento

ggplot()+
  geom_sf(data = zonas_censales)+
  geom_sf(data = recorridos)+
  geom_sf(data = recorridos_no_maximales)+
  geom_sf(data = paradas)+
  theme_bw()

Observamos que el STM tiene alcance mas alla de las fronteras del departamento, pero en este trabajo unicamente utilizaremos aquellas paradas que estan dentro del departamento, para esto hacemos un spatial join para poder filtrar las paradas que no estan en el departamento y tambien asignarles el codigo de la zona censal dentro de la que pertenecen, asi podemos empezar a trabajar con los datos de transporte publico dentro de cada zona censal.

Flitro de zonas y paradas a utilizar

sf_use_s2(FALSE)

## Spherical geometry (s2) switched off

paradas <- paradas %>%
  st_make_valid() %>%
  st_join(zonas_censales)

## although coordinates are longitude/latitude, st_intersects assumes that they are planar

## although coordinates are longitude/latitude, st_intersects assumes that they are planar

paradas <- paradas %>%
  filter(!is.na(CODSEG))

Observamos los datos de las secciones censales para ver en que nivel de detalle analizamos su subdivision

head(zonas_censales)

## Simple feature collection with 6 features and 13 fields
## Geometry type: MULTIPOLYGON
## Dimension:     XY
## Bounding box:  xmin: -56.20501 ymin: -34.91282 xmax: -56.12181 ymax: -34.89936
## Geodetic CRS:  WGS 84
##        AREA PERIMETER DEPTO SECCION SEGMENTO LOCALIDAD CODSEC CODSEG CODLOC
## 1  14476.50 1785.1361     1       0        0        20    100 100000   1020
## 2 122200.01 1563.4714     1       1        1        20    101 101001   1020
## 3 151608.67 1676.3479     1       1        2        20    101 101002   1020
## 4  99377.28 1456.9389     1       1        3        20    101 101003   1020
## 5  54395.51 1040.6595     1       1      104        20    101 101104   1020
## 6  41055.69  838.3599     1       1      105        20    101 101105   1020
##    NOMBDEPTO    NOMBLOC CDEPTO_ISO CLOC_ISO                       geometry
## 1 MONTEVIDEO MONTEVIDEO       UYMO  UYMOMON MULTIPOLYGON (((-56.13303 -...
## 2 MONTEVIDEO MONTEVIDEO       UYMO  UYMOMON MULTIPOLYGON (((-56.20449 -...
## 3 MONTEVIDEO MONTEVIDEO       UYMO  UYMOMON MULTIPOLYGON (((-56.1998 -3...
## 4 MONTEVIDEO MONTEVIDEO       UYMO  UYMOMON MULTIPOLYGON (((-56.19956 -...
## 5 MONTEVIDEO MONTEVIDEO       UYMO  UYMOMON MULTIPOLYGON (((-56.19834 -...
## 6 MONTEVIDEO MONTEVIDEO       UYMO  UYMOMON MULTIPOLYGON (((-56.19813 -...

zonas_censales <- zonas_censales %>%
  mutate(CODSEC = as.factor(CODSEC))

ggplot()+
  geom_sf(data = zonas_censales,aes(fill=CODSEC))+
  labs(title = "Secciones Departamentles",
       y="Latitud",
       x="Longitud",
       legend="Seccion",
       caption="Fuente: INE")+
  theme_bw()

Observamos que tiene 3 niveles de subdivision, el departamento, la seccion y el segmento. como el departamento solo usaremos uno observaremos el nivel de detalle de cada seccion para observar si con este nivel es suficiente para sacar conclusiones geograficas o es necesario ir un nivel mas abajo y analiza por segmento censal

Cramos un nuevo SF de poligonos por zonas Censales para facilitar el analisis

zonas_censales_CODSEC <- zonas_censales %>%
  group_by(CODSEC) %>%
  summarise()

## although coordinates are longitude/latitude, st_union assumes that they are planar
## although coordinates are longitude/latitude, st_union assumes that they are planar
## although coordinates are longitude/latitude, st_union assumes that they are planar
## although coordinates are longitude/latitude, st_union assumes that they are planar
## although coordinates are longitude/latitude, st_union assumes that they are planar
## although coordinates are longitude/latitude, st_union assumes that they are planar
## although coordinates are longitude/latitude, st_union assumes that they are planar
## although coordinates are longitude/latitude, st_union assumes that they are planar
## although coordinates are longitude/latitude, st_union assumes that they are planar
## although coordinates are longitude/latitude, st_union assumes that they are planar
## although coordinates are longitude/latitude, st_union assumes that they are planar
## although coordinates are longitude/latitude, st_union assumes that they are planar
## although coordinates are longitude/latitude, st_union assumes that they are planar
## although coordinates are longitude/latitude, st_union assumes that they are planar
## although coordinates are longitude/latitude, st_union assumes that they are planar
## although coordinates are longitude/latitude, st_union assumes that they are planar
## although coordinates are longitude/latitude, st_union assumes that they are planar
## although coordinates are longitude/latitude, st_union assumes that they are planar
## although coordinates are longitude/latitude, st_union assumes that they are planar
## although coordinates are longitude/latitude, st_union assumes that they are planar
## although coordinates are longitude/latitude, st_union assumes that they are planar
## although coordinates are longitude/latitude, st_union assumes that they are planar
## although coordinates are longitude/latitude, st_union assumes that they are planar
## although coordinates are longitude/latitude, st_union assumes that they are planar
## although coordinates are longitude/latitude, st_union assumes that they are planar
## although coordinates are longitude/latitude, st_union assumes that they are planar
## although coordinates are longitude/latitude, st_union assumes that they are planar

Alcance espacial

ggplot()+
  geom_sf(data = zonas_censales_CODSEC,aes(fill=CODSEC))+
  geom_rect(aes(xmin = -56.235, xmax = -56.11, ymin = -34.85, ymax = -34.94), color = "red", fill = NA)+
  labs(title = "Secciones Departamentles",
       y="Latitud",
       x="Longitud",
       legend="Seccion",
       caption="Fuente: INE")+
  theme_bw()

Para facilitar el calculo y el manejo de 22 millones de datos utilizaremos group by y summarise para resumir el dataset y poder graficar los datos por zona censal

datos_julio_resumido <- datos_julio %>%
  group_by(codigo_parada_origen,
           sevar_codigo,
           dsc_linea) %>%
  summarise(Total = sum(cantidad_pasajeros)) %>%
  rename(COD_UBIC_P = codigo_parada_origen) %>%
  mutate(parada_linea = paste(dsc_linea,sevar_codigo,COD_UBIC_P,sep="_")) %>%
  select(parada_linea,
         Total)

## `summarise()` has grouped output by 'codigo_parada_origen', 'sevar_codigo'. You
## can override using the `.groups` argument.
## Adding missing grouping variables: `COD_UBIC_P`, `sevar_codigo`

paradas <- paradas %>%
  mutate(parada_linea = paste(DESC_LINEA,COD_VARIAN,COD_UBIC_P,sep="_"))

unimos el dataset resumido al de paradas y zonas censales para poder graficarlos

paradas <- paradas %>%
  left_join(datos_julio_resumido,by = "parada_linea") %>%
  mutate_each(funs(replace(., which(is.na(.)), 0))) %>%
  mutate(CODSEC = as.character(CODSEC))

## Warning: `funs()` was deprecated in dplyr 0.8.0.
## Please use a list of either functions or lambdas: 
## 
##   # Simple named list: 
##   list(mean = mean, median = median)
## 
##   # Auto named with `tibble::lst()`: 
##   tibble::lst(mean, median)
## 
##   # Using lambdas
##   list(~ mean(., trim = .2), ~ median(., na.rm = TRUE))
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was generated.

## Warning: `mutate_each_()` was deprecated in dplyr 0.7.0.
## Please use `across()` instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was generated.

Graficamos el total de boletos vendidos por seccion

ggplot(paradas %>%
         st_set_geometry(NULL) %>%
         group_by(CODSEC) %>%
         summarise(Total = sum(Total)))+
  geom_bar(aes(y=reorder(CODSEC,(Total)),weight = Total),fill = "#60E6C8")+
  labs(title="Total de Viajes Vendidos por Seccion Departamental",
       y="Seccion",
       x="Total",
       caption="Fuente: STM")+
  theme_bw()

Observamos que las secciones 110 y 111 son las dos que mas viajes vendidos tienen, ahora observamos espacialmente

viajes_por_seccion <- paradas %>%
  st_set_geometry(NULL) %>%
  group_by(CODSEC) %>%
  summarise(Total = sum(Total))

zonas_censales_CODSEC <- left_join(zonas_censales_CODSEC,viajes_por_seccion,by = "CODSEC")

Mapa del Total por Seccion Censal

ggplot()+
  geom_sf(data = zonas_censales_CODSEC,aes(fill=Total))+
  scale_fill_gradient(low = "#FDFF62",high = "#FF6262")+
  labs(title = "Total de Viajes Vendidos por Zona Censal",
       x = "Longitud",
       y = "Latitud",
       legend = "Viajes Vendidos",
       caption = "Fuente: STM")+
  theme_bw()

Se puede observar que la zona sureste de la ciudad es la seccion con mas viajes vendidos, superando los dos millones de viajes en el mes de Julio 2022, para concluir aun mas sobre estas zonas se puede hacer mas hincapie en relacion a la poblacion y puestos de trabajo en cada zona censal, pero sobre esto se indagara en los proximos trabajos practicos.

Trabajo Practico N° 2

Elaborado por Ramiro Duarte y Gonzalo Machin

se puede acceder a este trabajo practico tambien en el siguiente Link

Alcance geografico del Sistema de Transporte Metropolitano

En la segunda entraga de este trabajo practico abordaremos que tan bueno es el alcance del transporte publico y entenderemos si todos los rincones del departamento tienen la misma accesibilidad al sistema

Para esto calcularemos la distancia a la parada mas cercana de diferentes puntos de las zonas censales

Primero responderemos una pregunta que quedo pendiente del TP anterior, existe una relacion entre la cantidad de viajes vendidos y la densidad de poblacion o la superficie?

Estudio de la Poblacion

Para esto necesitaremos los datos del censo 2011

Poblacion <- read.csv("data/censo/personas_por_zona.csv")

Con los datos de poblacion cargados, procedemos a limpiarlos un poco y hacer un join para asignarlos a cada seccion censal

Poblacion <- Poblacion %>%
  mutate(CODSEG = ((DPTO*100000) + (SECC*1000) + SEGM)) %>%
  mutate(CODSEC = ((DPTO*100) + SECC))

Poblacion_por_seccion <- Poblacion %>%
  group_by(CODSEC) %>%
  summarise(Total = sum(Total))

zonas_censales_CODSEC <- zonas_censales_CODSEC %>%
  mutate(CODSEC = as.numeric(CODSEC)) %>%
  left_join(Poblacion_por_seccion,by = "CODSEC")

zonas_censales_CODSEC <- zonas_censales_CODSEC %>%
  rename(Viajes_totales = Total.x,
         Poblacion = Total.y) 

zonas_censales_CODSEC <- zonas_censales_CODSEC %>%
  filter(!is.na(Poblacion))

Densidad de poblacion

Calculamos el area de cada zona censal y la densidad de poblacion

zonas_censales_CODSEC <- zonas_censales_CODSEC %>%
  mutate(Area_km2 = units::set_units(st_area(zonas_censales_CODSEC), km^2)) %>%
  mutate(Area_km2 = sub(" [km^2].*","",Area_km2)) %>%
  mutate(Area_km2 = as.numeric(Area_km2)) %>%
  mutate(Densidad = Poblacion/Area_km2)

pob <- ggplot(zonas_censales_CODSEC,aes(x=Poblacion,y=Viajes_totales))+
  geom_point()+
  geom_smooth(method=lm, se=FALSE, col='red', size=0.5)+
  labs(title = "Viajes Totales, Poblacion y Densidad",
       x="Poblacion",
       y = "Viajes Vendidos")+
  theme_bw()

density <- ggplot(zonas_censales_CODSEC,aes(x=Densidad,y=Viajes_totales))+
  geom_point()+
  geom_smooth(method=lm, se=FALSE, col='red', size=0.5)+
  labs(x="Densidad de Poblacion",
       y = "Viajes Vendidos")+
  theme_bw()

pob + density

## `geom_smooth()` using formula 'y ~ x'
## `geom_smooth()` using formula 'y ~ x'

Se puede observar una mayor correlacion y relacion lineal de la poblacion con la cantidad de viajes vendidos que la densidad.

Ahora estudiaremos la accesibilidad del Transporte Publico en la ciudad, para esto cargaremos primero las direcciones oficiales de la ciudad

Accesibilidad

Direcciones

direcciones <- st_read("data/direcciones/v_mdg_accesos.shp", options = "ENCODING=UTF-8")

## options:        ENCODING=UTF-8 
## Reading layer `v_mdg_accesos' from data source 
##   `E:\Gonzalo Machin\Maestria\Materias\CDC 2\tp1_ciencia_de_datos_2\data\direcciones\v_mdg_accesos.shp' 
##   using driver `ESRI Shapefile'
## Simple feature collection with 372443 features and 7 fields
## Geometry type: POINT
## Dimension:     XY
## Bounding box:  xmin: 552746.8 ymin: 6134530 xmax: 588662.6 ymax: 6159661
## Projected CRS: WGS 84 / UTM zone 21S

Observamos el dataset cargado

ggplot()+
  geom_sf(data = zonas_censales_CODSEC)+
  geom_sf(data = direcciones)+
  theme_bw()

Ahora volveremos a transformar a CRS 32721 el .shp de paradas

paradas <- paradas %>%
  st_transform(32721)

Para estudiar la accesibilidad crearemos un buffer de 500 metros de cada parada para saber que tantas lineas y variantes tiene acceso una direccion particular, en este caso cabe destacar que nuestro dataset de paradas crea un punto de iguales coordenadas para cada variante de linea de omnibus que para en esa parada, por lo tanto estariamos calculando el total de variantes a las que se tiene acceso desde una direccion particular..

ST_BUFFER

paradas_buffer <- st_buffer(paradas, 500)

ggplot()+
  geom_sf(data = zonas_censales_CODSEC)+
  geom_sf(data = paradas_buffer)+
  geom_sf(data = paradas)+
  theme_bw()

Con el buffer podemos saber cuantas direcciones son abastecidas y a cuantas lineas tienen acceso, dado que queremos conservar aquellas direcciones que no tienen acceso al STM a menos de 500 mts, utilizaremos st_join para unir ambos datasets

Primero uniremos los buffers para mejorar la velocidad de calculo del st_join

paradas_buffer_lineas <- paradas_buffer %>%
  group_by(DESC_LINEA) %>%
  summarise(Paradas = n_distinct(COD_UBIC_P.x))

ggplot()+
  geom_sf(data = zonas_censales_CODSEC)+
  geom_sf(data = paradas_buffer_lineas)+
  geom_sf(data = paradas)+
  theme_bw()

Union espacial y calculo de la accesibilidad

Accesibilidad <- st_join(direcciones,paradas_buffer_lineas)

Ahora calculamos el total de Lineas a las que una direccion particular puede acceder a menos de 500 mts

Accesibilidad <- Accesibilidad %>%
  mutate(Direccion = paste(NOM_CALLE,NUM_PUERTA,sep=" ")) %>%
  mutate_all(~replace(., is.na(.), 0)) %>%
  mutate(Accesos = if_else(Paradas == 0,0,1))

Direcciones_accesibles <- Accesibilidad %>%
  st_set_geometry(NULL) %>%
  group_by(Direccion) %>%
  summarise(Lineas_accesibles = sum(Accesos))

Ahora procedemos a unir y graficar el total de paradas y lineas accesibles a menos de 500 metros desde cada direccion de la ciudad

direcciones <- direcciones %>%
  mutate(Direccion = paste(NOM_CALLE,NUM_PUERTA,sep=" ")) %>%
  left_join(Direcciones_accesibles,by = "Direccion")

my_breaks = c(0,2,5,15,50,200,500)

ggplot()+
  geom_sf(data = zonas_censales_CODSEC)+
  geom_sf(data = paradas_buffer_lineas)+
  geom_sf(data = direcciones,aes(color = Lineas_accesibles))+
  labs(title = "Lineas y Paradas Accesibles desde cada Direccion",
       x = "Longitud",
       y = "Latitud",
       legend = "Accesibilidad")+
  scale_color_gradient(trans = "log",breaks = my_breaks,labels = my_breaks)+
  theme_bw()

## Warning: Transformation introduced infinite values in discrete y-axis

Se puede apreciar que aquellas direcciones que se ubican cerca o sobre de las principales arterias de la ciudad tienen una mayor accesibilidad al STM, mientras que en gris se pueden observar aquellos lugares que estan por fuera, predominantemente ubicadas en el norte y oeste de la ciudad, esto se puede explicar probablemente por la poblacion y baja densidad de poblacion.

Para observar mejor el impacto de las arterias de circulacion sobre la accesibilidad podemos hacer un hexbin map para visualizar mejor

Hexbin de accesibilidad

montevideo_grid <- st_make_grid(zonas_censales_CODSEC,
                             n=c(125,125),
                             what = 'polygons',
                             square = FALSE,
                             flat_topped = TRUE) %>%
  st_as_sf() %>%
  mutate(area = st_area(.)) %>%
  mutate(ID = row_number())

ggplot()+
  geom_sf(data = montevideo_grid)+
  geom_sf(data = zonas_censales_CODSEC,alpha = 0.4)+
  theme_bw()

Direcciones en cada hexbin y calculo del promedio de accesibilidad dentro de cada cuadrante

direcciones <- direcciones %>%
  st_transform(4326)

direcciones_hex <- st_join(direcciones,montevideo_grid)

## although coordinates are longitude/latitude, st_intersects assumes that they are planar
## although coordinates are longitude/latitude, st_intersects assumes that they are planar

hex_accesibilidad <- direcciones_hex %>%
  st_set_geometry(NULL) %>%
  group_by(ID) %>%
  summarise(Lineas_accesibles = floor(mean(Lineas_accesibles)))

Join y grafico de Hex

montevideo_grid <- left_join(montevideo_grid,hex_accesibilidad,by = "ID")

montevideo_grid <- montevideo_grid %>%
  filter(!is.na(Lineas_accesibles))

ggplot()+
  geom_sf(data = zonas_censales_CODSEC,fill = NA)+
  geom_sf(data = montevideo_grid,aes(fill=Lineas_accesibles))+
  scale_fill_gradientn(colours = c("#7FF2A4","#CAF087","#F6F296","#F2B637","#CD4206","#5F1F03"),
                       guide = guide_colourbar(barheight = 10))+
  labs(title = "Accesibilidad al STM por sector de la Ciudad",
       x = "Longitud",
       y = "Latitud",
       legend = "Accesibilidad")+
  theme_bw()

Trabajo Practico N° 3

Elaborado por Ramiro Duarte y Gonzalo Machin

se puede acceder a este trabajo practico tambien en el siguiente Link

Representación de mapas con OpenStreetMap

Descarga de mapa base de la Ciudad de Montevideo.

Descargamos el cuadro delimitador mediante la función getbb.

bbox_montevideo <- getbb("Montevideo, Uruguay")

Mostramos el resultado.

bbox_montevideo

##         min       max
## x -56.43140 -56.02250
## y -34.93806 -34.70185

Descargamos el mapa base mediante la función get_stamenmap.

mapa_montevideo <- get_stamenmap(bbox=bbox_montevideo,
              maptype="toner-lite",
              zoom=12)

## Source : http://tile.stamen.com/toner-lite/12/1405/2469.png

## Source : http://tile.stamen.com/toner-lite/12/1406/2469.png

## Source : http://tile.stamen.com/toner-lite/12/1407/2469.png

## Source : http://tile.stamen.com/toner-lite/12/1408/2469.png

## Source : http://tile.stamen.com/toner-lite/12/1409/2469.png

## Source : http://tile.stamen.com/toner-lite/12/1410/2469.png

## Source : http://tile.stamen.com/toner-lite/12/1405/2470.png

## Source : http://tile.stamen.com/toner-lite/12/1406/2470.png

## Source : http://tile.stamen.com/toner-lite/12/1407/2470.png

## Source : http://tile.stamen.com/toner-lite/12/1408/2470.png

## Source : http://tile.stamen.com/toner-lite/12/1409/2470.png

## Source : http://tile.stamen.com/toner-lite/12/1410/2470.png

## Source : http://tile.stamen.com/toner-lite/12/1405/2471.png

## Source : http://tile.stamen.com/toner-lite/12/1406/2471.png

## Source : http://tile.stamen.com/toner-lite/12/1407/2471.png

## Source : http://tile.stamen.com/toner-lite/12/1408/2471.png

## Source : http://tile.stamen.com/toner-lite/12/1409/2471.png

## Source : http://tile.stamen.com/toner-lite/12/1410/2471.png

## Source : http://tile.stamen.com/toner-lite/12/1405/2472.png

## Source : http://tile.stamen.com/toner-lite/12/1406/2472.png

## Source : http://tile.stamen.com/toner-lite/12/1407/2472.png

## Source : http://tile.stamen.com/toner-lite/12/1408/2472.png

## Source : http://tile.stamen.com/toner-lite/12/1409/2472.png

## Source : http://tile.stamen.com/toner-lite/12/1410/2472.png

Mostramos el mapa con la función ggmap.

ggmap(mapa_montevideo)

Delimitamos el polígono de la ciudad.

poligono_montevideo <- getbb("Montevideo, Uruguay",
                             format_out="sf_polygon")

Vemos el resultado.

poligono_montevideo

## Simple feature collection with 2 features and 0 fields
## Geometry type: POLYGON
## Dimension:     XY
## Bounding box:  xmin: -56.4314 ymin: -34.93806 xmax: -56.0225 ymax: -34.70185
## Geodetic CRS:  WGS 84
##                         geometry
## 1 POLYGON ((-56.4314 -34.8281...
## 2 POLYGON ((-56.4314 -34.8281...

ggmap(mapa_montevideo)+
  geom_sf(data=poligono_montevideo,
          fill=NA,
          color="darkgreen",
          size=1,
          inherit.aes=FALSE)+
  labs(title="Montevideo",
       subtitle="Uruguay",
       caption="Fuente: Open Street Map")+
  theme_void()

## Coordinate system already present. Adding new coordinate system, which will replace the existing one.

Descarga “líneas”

Para mantener coherencia con el resto de los trabajos prácticos de la materia, descargamos las calles de la ciudad.

En primer lugar, asignamos el boundig box.

calles_montevideo <- opq(bbox_montevideo) %>%
  add_osm_feature(key = "highway",
                  value = c("motorway",
                          "trunk",
                          "primary",
                          "secondary",
                          "tertiary",
                          "unclassified",
                          "residential"))

Utilizamos la función osmdata_sf para descargar la informació y la mostramos.

calles_montevideo <- osmdata_sf(calles_montevideo)

Nos quedamos con las líneas.

calles_montevideo <- calles_montevideo$osm_lines

Graficamos el resultado.

ggplot()+
  geom_sf(data=poligono_montevideo, fill=NA, color="darkgreen", size=1)+
    geom_sf(data=calles_montevideo)+
  labs(title="Calles",
       subtitle="Montevideo, Uruguay",
       caption="Fuente: Open Street Map")+
  theme_void()

Lo unimos con nuestro mapa base.

ggmap(mapa_montevideo)+
  geom_sf(data=poligono_montevideo, fill=NA, color="darkgreen", size=1, inherit.aes = FALSE)+
    geom_sf(data=calles_montevideo, inherit.aes = FALSE)+
  labs(title="Calles",
       subtitle="Montevideo, Uruguay",
       caption="Fuente: Open Street Map")+
  theme_void()

## Coordinate system already present. Adding new coordinate system, which will replace the existing one.

dim(calles_montevideo)

## [1] 14733   114

calles_montevideo %>%
  group_by(highway) %>%
  summarise(cantidad=n())

## although coordinates are longitude/latitude, st_union assumes that they are planar
## although coordinates are longitude/latitude, st_union assumes that they are planar
## although coordinates are longitude/latitude, st_union assumes that they are planar
## although coordinates are longitude/latitude, st_union assumes that they are planar
## although coordinates are longitude/latitude, st_union assumes that they are planar
## although coordinates are longitude/latitude, st_union assumes that they are planar

## Simple feature collection with 6 features and 2 fields
## Geometry type: MULTILINESTRING
## Dimension:     XY
## Bounding box:  xmin: -56.44871 ymin: -34.93419 xmax: -55.99676 ymax: -34.65824
## Geodetic CRS:  WGS 84
## # A tibble: 6 x 3
##   highway      cantidad                                                 geometry
##   <chr>           <int>                                    <MULTILINESTRING [°]>
## 1 primary           620 ((-56.08889 -34.73053, -56.08974 -34.7317, -56.09039 -3~
## 2 residential     11144 ((-56.03465 -34.84728, -56.03397 -34.84764, -56.03344 -~
## 3 secondary         842 ((-56.11317 -34.72069, -56.11315 -34.72177, -56.11308 -~
## 4 tertiary         1458 ((-56.22758 -34.72739, -56.22695 -34.72737, -56.22586 -~
## 5 trunk             447 ((-56.35509 -34.78715, -56.35661 -34.78021), (-56.10173~
## 6 unclassified      222 ((-56.17819 -34.71774, -56.17826 -34.71765, -56.1795 -3~

Hay 14.733 observaciones, dividiendose en: - 620 primarias. - 11.144 residenciales. - 842 secundarias. - 1.458 terciarias. - 447 troncales. - 222 sin clasificar.

Eliminamos los tramos fuera de nuestra área de influencia.

calles_montevideo <- st_intersection(calles_montevideo, poligono_montevideo)

## although coordinates are longitude/latitude, st_intersection assumes that they are planar

## Warning: attribute variables are assumed to be spatially constant throughout all
## geometries

dim(calles_montevideo)

## [1] 24056   114

calles_montevideo %>%
  group_by(highway) %>%
  summarise(cantidad=n())

## although coordinates are longitude/latitude, st_union assumes that they are planar
## although coordinates are longitude/latitude, st_union assumes that they are planar
## although coordinates are longitude/latitude, st_union assumes that they are planar
## although coordinates are longitude/latitude, st_union assumes that they are planar
## although coordinates are longitude/latitude, st_union assumes that they are planar
## although coordinates are longitude/latitude, st_union assumes that they are planar

## Simple feature collection with 6 features and 2 fields
## Geometry type: MULTILINESTRING
## Dimension:     XY
## Bounding box:  xmin: -56.43103 ymin: -34.93419 xmax: -56.03084 ymax: -34.70411
## Geodetic CRS:  WGS 84
## # A tibble: 6 x 3
##   highway      cantidad                                                 geometry
##   <chr>           <int>                                    <MULTILINESTRING [°]>
## 1 primary          1126 ((-56.19722 -34.91148, -56.19859 -34.9111), (-56.19859 ~
## 2 residential     17918 ((-56.19148 -34.72162, -56.1917 -34.72148), (-56.18441 ~
## 3 secondary        1640 ((-56.06011 -34.86677, -56.05993 -34.86667), (-56.05993~
## 4 tertiary         2488 ((-56.21423 -34.74944, -56.21422 -34.74923), (-56.21422~
## 5 trunk             692 ((-56.35509 -34.78715, -56.35577 -34.78405), (-56.10173~
## 6 unclassified      192 ((-56.28868 -34.72325, -56.2886 -34.72319), (-56.2886 -~

Hay 14.733 observaciones, dividiendose en: - X primarias. - X residenciales. - X secundarias. - X terciarias. - X troncales. - X sin clasificar.

Eliminamos los tramos fuera de nuestra área de influencia.

Graficamos.

ggmap(mapa_montevideo)+
  geom_sf(data=poligono_montevideo, fill=NA, color="darkgreen", size=1, inherit.aes = FALSE)+
    geom_sf(data=calles_montevideo, inherit.aes = FALSE)+
  labs(title="Calles",
       subtitle="Montevideo, Uruguay",
       caption="Fuente: Open Street Map")+
  theme_void()

## Coordinate system already present. Adding new coordinate system, which will replace the existing one.

Descarga “puntos”

Descargamos las paradas de omnibus.

En primer lugar, asignamos el boundig box.

paradas_montevideo <- opq(bbox_montevideo) %>%
  add_osm_feature(key = "highway",
                  value = "bus_stop")

Utilizamos la función osmdata_sf para descargar la información

paradas_montevideo <- osmdata_sf(paradas_montevideo)

paradas_montevideo

## Object of class 'osmdata' with:
##                  $bbox : -34.938056,-56.4313997,-34.7018526,-56.0225006
##         $overpass_call : The call submitted to the overpass API
##                  $meta : metadata including timestamp and version numbers
##            $osm_points : 'sf' Simple Features Collection with 204 points
##             $osm_lines : NULL
##          $osm_polygons : 'sf' Simple Features Collection with 1 polygons
##        $osm_multilines : NULL
##     $osm_multipolygons : NULL

Nos quedamos con los valores “punto”.

paradas_montevideo <- paradas_montevideo$osm_points

Existen 204 observaciones.

Graficamos.

ggmap(mapa_montevideo)+
  geom_sf(data=poligono_montevideo, fill=NA, color="darkgreen", size=1, inherit.aes = FALSE)+
    geom_sf(data=paradas_montevideo, inherit.aes = FALSE, aes(color=highway))+
  labs(title="Paradas",
       subtitle="Montevideo, Uruguay",
       color="Tipo",
       caption="Fuente: Open Street Map")+
  scale_color_manual(values="darkorange")+
  theme_void()

## Coordinate system already present. Adding new coordinate system, which will replace the existing one.

Se ve que hay puntos por fuera del área de estudio. Se verá primero la dimensión y luego se los quitará con la función st_intersection.

dim(paradas_montevideo)

## [1] 204  35

paradas_montevideo %>%
  group_by(highway) %>%
  summarise(cantidad=n())

## although coordinates are longitude/latitude, st_union assumes that they are planar
## although coordinates are longitude/latitude, st_union assumes that they are planar

## Simple feature collection with 2 features and 2 fields
## Geometry type: MULTIPOINT
## Dimension:     XY
## Bounding box:  xmin: -56.42333 ymin: -34.92366 xmax: -56.02256 ymax: -34.70244
## Geodetic CRS:  WGS 84
## # A tibble: 2 x 3
##   highway  cantidad                                                     geometry
##   <chr>       <int>                                             <MULTIPOINT [°]>
## 1 bus_stop      200 ((-56.42333 -34.75531), (-56.42215 -34.75554), (-56.40782 -~
## 2 <NA>            4 ((-56.22367 -34.72801), (-56.22366 -34.72779), (-56.22364 -~

Vemos que hay 200 paradas y 4 N/A.

paradas_montevideo <- st_intersection(paradas_montevideo, poligono_montevideo)

## although coordinates are longitude/latitude, st_intersection assumes that they are planar

## Warning: attribute variables are assumed to be spatially constant throughout all
## geometries

ggmap(mapa_montevideo)+
  geom_sf(data=poligono_montevideo, fill=NA, color="darkgreen", size=1, inherit.aes = FALSE)+
    geom_sf(data=paradas_montevideo, inherit.aes = FALSE, aes(color=highway))+
  labs(title="Paradas",
       subtitle="Montevideo, Uruguay",
       color="Tipo",
       caption="Fuente: Open Street Map")+
  scale_color_manual(values="darkorange")+
  theme_void()

## Coordinate system already present. Adding new coordinate system, which will replace the existing one.

paradas_montevideo %>%
  group_by(highway) %>%
  summarise(cantidad=n())

## although coordinates are longitude/latitude, st_union assumes that they are planar

## Simple feature collection with 1 feature and 2 fields
## Geometry type: MULTIPOINT
## Dimension:     XY
## Bounding box:  xmin: -56.35556 ymin: -34.92366 xmax: -56.03415 ymax: -34.76106
## Geodetic CRS:  WGS 84
## # A tibble: 1 x 3
##   highway  cantidad                                                     geometry
##   <chr>       <int>                                             <MULTIPOINT [°]>
## 1 bus_stop      204 ((-56.35556 -34.82972), (-56.35414 -34.7878), (-56.35386 -3~